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ABSTRACT 

We study the formation and evolution of Hii regions around the first stars formed at redshifts 
z = 10 — 30. We use a one-dimensional Lagrangian hydrodynamics code which self-consistently 
incorporates radiative transfer and non-equilibrium primordial gas chemistry. The star-forming region 
is defined as a spherical dense molecular gas cloud with a Population III star embedded at the center. 
We explore a large parameter space by considering, as plausible early star-forming sites, dark matter 
halos of mass Afhaio — 10^ — 10^ Mq, gas density profiles with a power-law index w = 1.5 — 2.25, and 
metal- free stars of mass Mgtar = 25 — 500 Mq . The formation of the Hii region is characterized by 
initial slow expansion of a weak D-type ionization front near the center, followed by rapid propagation 
of an R-type front throughout the outer gas envelope. We find that the transition between the two 
front types is indeed a critical condition for the complete ionization of halos of cosmological interest. 
In small mass (< 10^ M©) halos, the transition takes place within a few 10^ yr, yielding high escape 
fractions (> 80%) of both ionizing and photodissociating photons. The gas is effectively evacuated 
by a supersonic shock, with the mean density within the halo decreasing to < Icm^'^ in a few million 
years. In larger mass (> 10^ M©) halos, the ionization front remains to be of D-type over the lifetime 
of the massive star, the Hii region is confined well inside the virial radius, and the escape fractions are 
essentially zero. We derive an analytic formula, that reproduces well the results of our simulations, for 
the critical halo mass below which the gas is completely ionized. We discuss immediate implications 
of the present results for the star formation history and early reionization of the Universe. 
Subject headings: cosmology: theory - reionization of the Universe - stars: Population III - radiative 
transfer 



1. INTRODUCTION 

The emergence of the first generation stars have signif- 
icant impacts on the thermal state of the inter-galactic 
medium (IGM) in the early universe. The initially neu- 
tral cosmic gas is photoionized and photoheated by radi- 
ation from the first stars. This so-called radiative feed- 
back from the first stars is of considerable cosmological 
interest. It can not only self-regulate the first star for- 
mation, but also affect the formation and evolution of 
proto-galaxies. 

Theories based on Cold Dark Matter (CDM) predict 
that the first cosmological objects form when small mass 
(~ 10^ Mq) dark halos assemble at redshifts z = 20 — 30 
(Couchman & Rees 1986; Tegmark et al. 1997; Abel et 
al. 1998; Yoshida et al. 2003a). Stars formed in these 
pre-galactic objects may have substantially contributed 
to cosmic reionization. Recent measurement of the large 
Thomson optical depth by the WMAP satellite indicates 
that the universe was reionized at an epoch as early as 
z ~ 17 (Kogut et al. 2003; Spergel et al. 2003), sup- 
porting the above notion that the first generation stars 



formed at z > 20. Theoretical studies of the processes of 
cosmic reionization by stellar sources suggest that reion- 
ization proceeds first by the formation of individual Hii 
regions around radiation sources (galaxies) , and then by 
percolation of the growing Hii bubbles (Gnedin & Os- 
triker 1997; Ricotti et al. 2002; Sokasian et al. 2004). 
The shape and the extension of the individual Hii regions 
critically determine the global topology of the ionized re- 
gions in a cosmological volume at different epochs during 
reionization. 

Studies on the formation of Hii regions in dense gas 
clouds date back to the seminal work by Stromgren 
(1939). Since then the structure of Hii regions and the 
interaction with the surrounding medium have been ex- 
tensively studied (see Yorke 1986 for a review). An im- 
portant advance has been made by Franco, Tenorio-Tagle 
& Bodenheimer (1990) who considered realistic condi- 
tions in which the initial gas density profile is given by 
a power-law. Franco et al. showed that a hydrodynamic 
shock effectively sweeps the ambient medium into a thin 
shell and the gas density profile is significantly modified 
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from the initial power-law shape. Shu et al. (2002) ob- 
tained the self-similar solutions for the expansion of self- 
gravitating Hii regions. On cosmological backgrounds, 
Shapiro & Giroux (1987) studied the evolution of the 
cosmological Stromgren sphere around luminous quasars 
in an expanding universe. Ricotti & ShuU (2000) com- 
puted the UV photon consumption in small mass halos to 
estimate the photon escape fraction. Their calculations, 
however, do not include the hydrodynamic response of 
the photoheated gas and are not adequate to address 
the dynamical evolution of Hii regions. More recently, 
Whalen, Abel & Norman (2004) carried out a numer- 
ical simulation of ionization front (I-front) propagation 
starting from a realistic initial density profile for the first 
star-forming cloud. Whalen et al. found that the escape 
fraction of ionizing photons is close to unity in a partic- 
ular case with halo mass 7 x 10^ M© at z ~ 20. It is by 
no means trivial if such high escape fraction is achieved 
for different masses and redshifts. 

In the present paper, we study the evolution of Hii 
regions around Population III stars for a wide range of 
halo mass, redshift, and gas density profile. Specifically, 
we consider an "inside-out" situation, where the central 
source ionizes the surrounding gas and drives a super- 
sonic gas flow outward. While the overall configuration 
of our simulations may appear similar to those in con- 
ventional studies on Hii regions, we have several critical 
ingredients relevant to the case of early generation star- 
formation; 1) gravitational force exerted by the host dark 
matter halo, 2) radiation force by the central star, 3) ra- 
diative transfer of UV photons, 4) chemical reactions in- 
cluding formation and destruction of hydrogen molecules, 
and 5) radiative cooling and heating processes, includ- 
ing inverse-Compton cooling particularly important at 
z> 10. 

Our inside-out simulations are complementary to a 
number of works that focused on radiative feedback from 
an external field (Umemura & Ikeuchi 1984; Rees 1986; 
Bond, Szalay & Silk 1988; Efstathiou 1992; Thoul & 
Weinberg 1996; Kepner, Babul & Spergel; Barkana & 
Loeb 1999; Kitayama & Ikeuchi 2000; Susa & Umemura 
2000, 2004; Kitayama et al. 2000, 2001; Shapiro, Iliev 
& Raga 2003, 2004). While external radiation can easily 
ionize and blow away the outer envelope of gas clouds, 
photo-evaporation proceeds less effectively in the dens- 
est central part. Shapiro et al. (2004) indeed show 
that evaporation of an initially hydrostatic mini-halo 
(~ 10^ M0) takes about 100 million years when irradi- 
ated by a luminous source at a distance of ~ 1 Mpc. 
At the very center of the halo, self-shielding and dynam- 
ical infall of the gas can prevent complete evaporation 
and promote star formation (Kitayama et al. 2001; Susa 
& Umemura 2004). The subsequent evolution of the gas 
cloud will then be largely regulated by the feedback from 
the central star. As will be shown in the present paper, 
ionization of the ambient gas by an internal source is 
more rapid than that by the external radiation because 
of the greater incident radiation flux and the presence of 
the gas density gradient. 

The final state of the ionized gas near the central region 
is of particular interest in the study of the first super- 
nova explosions (e.g. Bromm, Yoshida & Hernquist 2003; 
Wada & Venkatesan 2003). The subsequent evolution of 
the ionized gas has also an important implication for the 



formation of proto-galaxies (Oh & Haiman 2003). For 
the gas clouds hosted by mini- halos with mass ~ 10^ Mq, 
one may naively expect that the gas is completely pho- 
toionized and eventually photoevaporated, because the 
virial temperature of the system is much lower than the 
characteristic temperature of photoionized gas. On the 
other hand, larger dark matter halos generate a deeper 
gravitational potential well, and thus it can effectively 
trap the hot, ionized gas within a small radius, rather 
than letting it move outward by pressure. The actual 
situation is likely to be far more complicated, depend- 
ing upon the luminosity and the spectrum of the radi- 
ation source, and also on the initial gas density profile. 
It is clearly important to make reliable predictions, un- 
der various physical conditions, as to how the gas is re- 
distributed by radiation from the central star. 

Throughout the present paper, we work with a A- 
dominated CDM cosmology with the matter density 
VLm = 0.3, the cosmological constant Qa = 0.7, the Hub- 
ble constant h = 0.7, and the baryon density VL-q = 0.05. 

2. THE SIMULATIONS 
2.1. Numerical Scheme 

We study the evolution of Hii regions around a primor- 
dial star using the radiation-hydrodynamics code of Ki- 
tayama etal. (2001). The code employs the second-order 
Lagrangian finite-difference scheme in spherically sym- 
metric geometry (Bowers & Wilson 1991; see also Thoul 
& Weinberg 1995). It treats self-consistently gravita- 
tional force, hydrodynamics, non-equilibrium chemistry 
of primordial gas, and the radiative transfer of UV pho- 
tons. In the present paper, we further incorporate the 
radiation force. We also adopt an artificial viscosity for- 
mulation of Caramana, Shashkov & Whalen (1998), de- 
signed to distinguish between shock-wave and uniform 
compression using an advection limiter. 

The basic equations are given by 



(1) 

/rad, (2) 

(3) 
(4) 



where r, m, p, p, u, and Aftot(< t) are the radius, mass, 
density, pressure, internal energy per unit mass, and the 
total mass inside r, respectively. The radiative heating 
and cooling rates per unit volume, Ti. and C, and the 
radiation force per unit mass, /rad, depend on the solu- 
tions of the radiative transfer equation and the chemical 
reactions. We therefore take the following procedure at 
each timestep we advance with the momentum equation 
(eq. 0). 

First, the direction and frequency dependent radiative 
transfer is worked out as described in detail in Appendix 
A. We solve both absorption and emission of ionizing 
(> 13.6 eV) photons and take account of self-shielding 
of H2 against the Lyman- Werner (LW, 11.2-13.6 eV) 
band photons following Draine & Bertoldi (1996). This 
yields the UV heating rate H and the radiation force /rad , 
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TABLE 1 

Properties of metal-free stars (Schaerer 2002). 



A/,tar [Mq] tiifc [Myr] T,ff [K] 7Vio„ [s"^] 



500 


1.90 


1.07 X 10^ 6.80 X 10^° 


200 


2.20 


9.98 X lO* 2.62 x lO^" 


80 


3.01 


9.33 X lO* 7.73 x 10*9 


25 


6.46 


7.08 X lO"* 7.58 X 10"'* 



together with the coefficients for photoioniziation of H, 
photodissociation of H2 and H^, and photodetachment 
of H~ . The obtained coefficients are used in the chemical 
reaction network mentioned below. 

Secondly, non-equilibrium chemical reactions are 
solved by an implicit scheme developed by Susa & Ki- 
tayama (2000) for the species e, H, H"*", H~, H2, and 
H^. ^ Unless otherwise stated, the reactions and the 
rates are taken from the compilation of Galli & Palla 
(1998). 

Finally, the energy equation (eq. [3|) is solved includ- 
ing UV heating and radiative cooling due to coUisional 
ionization, coUisional excitation, recombination, thermal 
bremsstrahlung, Compton scattering with the cosmic mi- 
crowave background (CMB) radiation, and rotational- 
vibrational excitation of H2. The atomic cooling rates 
are taken from the compilation of Fukugita & Kawasaki 
(1994) and molecular cooling rates from Galli & Palla 
(1998). 

The above procedure is repeated iteratively until the 
abundance of each species and the internal energy in each 
mesh converge within an accuracy of 0.1%. We have val- 
idated the accuracy of our code by carrying out a suite of 
conventional tests. The results of a "Stromgren sphere" 
test are presented in Appendix B. 

2.2. Central Source 

At the center of a gas cloud, we place a metal-free Pop- 
ulation III star. Its parameters are taken from Schaerer 
(2002) and listed in Tabled We vary the stellar mass, 
denoted by Mgtar, over the range 25 ~ 500 M0 and ap- 
proximate the spectrum by a black-body with the effec- 
tive temperature Toff. The luminosity is normalized so 
that the number of ionizing photons emitted per second, 
iVjon, matches the time averaged photon flux given in Ta- 
ble 5 of Schaerer (2002). For simplicity, we neglect its 
spectral evolution. Unless stated otherwise, the simula- 
tions are performed during the main sequence lifetime, 
iiifc, of the central star. 

2.3. Initial Conditions 

We assume that a star-forming gas cloud is embed- 
ded in a dark matter halo with the NFW density profile 
(Navarro, Frenk & White 1997): 

1 



Pdm fx 



x{i + xy 



(5) 



^ In the present paper, we consider only the hydrogen compo- 
nent of the IGM. Including helium will raise the temperature of 
photoionized gas and may affect the over all evolution of Hll re- 
gions. This point will be discussed in H3.3I 



where x — r /r^ is the radius normalized by the scale 
radius T-g. We follow Bullock et al. (2001) to determine rs 
for a halo with total mass Mhaio collapsing at redshift Zc, 
by extrapolating their formula to the lower halo masses 
and the higher redshifts. The dark matter density is 
normalized so that the dark matter mass within the virial 
radius is equal to Mdm = (1— ^B/^M)-A^haio- We assume 
that the dark matter density profile remains unchanged, 
since the halo dynamical timescale is much longer than 
the lifetime of the massive star. 

For the gas in the halo, we adopt a power-law density 
profile 

pocr-"", (6) 

and vary the index w from 1.5 to 2.25. This power- 
law profile are motivated by the results of recent three 
dimensional simulations of the primordial gas cloud for- 
mation (Abel et al. 2002; Yoshida et al. 2003a) which 
show that the gas density profile around the first star- 
forming regions is well described with w '^ 2 over a wide 
range of radii. The gas density is normalized so that the 
total gas mass within the halo virial radius is equal to 
Mgas = (riB/f^M)-Mhaio — -Mstar- We determine the radius 
of the inner-most gas shell rjn such that 



'^h(< r-m) = lO'^cm ^ 



(7) 



where nH(< f) is the average hydrogen number density 
within r. This density is comparable to that of primor- 
dial molecular cloud cores or fragments found in numer- 
ical simulations (Abel et al. 2002; Bromm et al. 2002; 
Nakamura & Umemura 2002), and gives r^^ sufficiently 
smaller than the characteristic size of the simulated ha- 
los. Our primary interest hence lies in the propagation 
of UV photons after they escape out of the dense molec- 
ular cloud. ^ If a gas shell falls below 0.5rin, it is moved 
to the center and ignored in the rest of the simulation, 
except in the calculation of the gravitational force. In or- 
der to trace the propagation of an I- front for a sufficiently 
long time, the outer boundary is taken at 1 ~ 3 times 
the virial radius depending on the run. The gas shells 
are initially spaced such that the shell mass increases by 
a constant ratio, typically '^ 1%, between the adjacent 
shells. We have checked that a total of 300 radial bins 
from the center to the outermost shell is sufficient to pro- 
duce converged results. Since the star-forming clouds are 
likely to be in the course of collapse, we assign the initial 
velocity 



Winit('') = - 



2GA/tot(< r) 



(8) 



We assume that the gas is initially isothermal with 
temperature Tjnit, given by the minimum of the virial 
temperature Tvir and a reference temperature Tmin- We 
take Tmin — 500 K, corresponding to the gas which has 
cooled by molecular hydrogen cooling. The final results 

^ The infalling envelope may prevent the escape of UV photons 
from dense proto-stellar clouds. Omukai & Inutsuka (2001) showed 
that, in spherically symmetric cases, Hll regions around a massive 
star cannot expand under significant mass accretion. Mass accre- 
tion onto a dense proto-stellar cloud, however, could occur along 
an aspherical disk due to its angular momentum. For disk-like 
geometries, escape of the UV photons from the proximity of the 
massive star will be greatly enhanced than in the spherically sym- 
metric case. In the present paper, we focus on such cases that the 
UV photons could successfully escape out of a proto-stellar disk 
and propagate through the surrounding halo. 
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Fig. 1. — Structure of an Hll region around a massive star with A/star = 200 Mq inside a halo with Mhaio = 10® Mq and w = 2.0 at 
Zc = 20. Radial profiles are shown at t = (dotted lines), 1.8 X 10^ yr (dashed lines), and 2.2 X 10® yr (solid lines) for (a) hydrogen density, 
(b) temperature, (c) ratio of radiation force to gravitational force, (d) electron fraction, (e) H2 fraction, and (f) radial velocity. 



are found to be insensitive to the choice of Tn- 



be- 



cause, as long as Tmin ^ lO^'K, the gas is initially al- 
most neutral and the opacity to ionizing photons is very 
large. The initial abundances are taken to be consistent 
with the above choice of the initial temperature based 
on the hydrodynamical simulations of Kitayama et al. 
(2001); Xe = 10-^ Xh- = 10^1°, Xh2 = 10-3, and 
Xfj2+ = 10^^^, where X^ = ni/nn is the fraction of the 
species i with respect to the total number of hydrogen 
atoms. 

3. RESULTS 

3.1. I-front propagation 

We first describe in detail the results of our fiducial 
runs with A'/gtar = 200 M©, w = 2.0, and Zc = 20. Simu- 



lations with other sets of parameters are also carried out 
and the results are presented in due course. 

Figure n shows the structure of the Hii region at i = 
(initial), 1.8 x 10^, and 2.2 x 10^ yr for a low mass case 
with Mhaio = 10^ Mq. The evolution in the early stage 
is characterized by the propagation of a weak D-type 
front into the surroundings; shock precedes the I-front. 
Because of the steep density gradient, the I-front even- 
tually overtakes the shock and changes into R-type. The 
transition takes place at i ~ 3 x 10^ yr (see also Figure 
Eland discussions below). In this so-called "champagne" 
phase, the low-density envelope is promptly ionized (Wel- 
ter 1980) and the gas temperature rises to a few times 
10^ K. As the shock propagates at the speed ~ 30 km 
s^^, it reaches ^ 70 pc within the lifetime of the massive 
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star, iiifo=2.2 Myr. 

The spatial resolution of the present simulation is suffi- 
cient to resolve the cooling layer behind the shocks; e.g., 
the I-front at t = 1.8 x 10^ yr in Figure H is resolved 
with more than 10 shells. We have also performed the 
runs with halving and doubling the gas shell numbers 
and confirmed that the results are almost identical. 

We notice that there appears a thin shell of H2 just in 
front of the Hii region (Fig. ^) as pointed out by Ri- 
cotti et al. (2001). This is due to a positive feedback on 
the H2 formation by an enhanced electron fraction at the 
temperature below 10*K. As the Hii region expands, the 
temperature in the shell exceeds lO'^K and H2 molecules 
are dissociated by collisions with ions. There still appears 
a new H2 shell in front of the I-front, but in a different 
position (both Eulerian and Lagrangian) from the pre- 



vious one. The H2 shell is thus likely to be short-lived 
as a result of I-front propagation. There may also be a 
case that the H2 formation is promoted after the central 
star fades away and the gas cools to temperatures be- 
low 10'*K. Such a possibility will be explored later in this 
subsection. 

Interestingly, the radiation force dominates over the 
gravitational force inside the Hii region when the I-front 
radius is smaller than ~ 30 pc (Fig. Q^). This value is 
a factor of ~ 4 smaller than that indicated by Haehnelt 
(1995), and easily accounted for by the fact that there is 
about 3-5 times more mass contributed by dark matter 
than baryons near the center in our simulations. 

For a higher mass halo with 10^ M0 illustrated in Fig- 
ure 121 the outward gas motion near the center is even- 
tually reverted to an inflow. Correspondingly, the shock 
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Fig. 3. — Evolution of (a) shock front radius (thick lines) and I-front radius (thin lines, whenever separated from the shock), (b) central 
hydrogen density, (c) mean HI fraction within rvir, and (d) mean H2 fraction within rvir, for halos with Afjjaio = lO'' M© (solid lines) and 
10^ Mq (dashed lines) shown in Figures and |5] respectively. Horizontal dotted lines denote virial radii of halos with Mji^io = 10® Mq 
and lO'^ Mq, while the vertical one indicates the lifetime of the central star with A/star = 200 Mq. 



radius starts decreasing. This reversion is explained by a 
combination of the deep gravitational potential well near 
the halo center and the infall of the envelope gas. Notice 
that the ionization fraction decreases sharply at r ~ 1 
pc and t = 1.8 X 10^ yr. The infalling material piles 
on the shock front, enhancing recombination (Fig. EJi). 
This is however just in a transition phase; the piled gas 
rapidly falls back to the center. Thereafter, the Hii re- 
gion is kept trapped within a few parsec radius because 
the gas density remains high (> 10'*cm~^) and recombi- 
nation balances photoionization. We have carried out a 
simulation for the same set of parameters but assuming 
that the gas is initially static. Even in this case, the D- 
type front maintains over iiifc=2.2 Myr and reaches only 



^ 40 pc, still well inside the virial radius. 

Figure|31further compares the time evolution in the two 
cases presented in Figures ^ and 13 For definiteness, we 
take the shock front radius at the peak of gas density and 
the I-front radius at the position with Xhi = 0.1. The 
mass weighted mean fractions of HI and H2 are computed 
within the virial radius rvir- For Mhaio = 10^ M©, the I- 
front detaches from the shock at t ^ 3 x 10^ yr and 
then propagates rapidly, exceeding the halo virial radius. 
The halo gas is almost completely ionized in the first 
few hundred thousand years and then a large fraction of 
ionizing and the LW photons can escape from the halo. 
For Mhaio = 10'' Mq, on the other hand, the I-front stalls 
at t ^ 3 X 10^ yr and then moves back inward. A large 
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Fig. 4. — Escape fractions of ionizing photons (> 13.6 eV, top 
panel) and the LW photons (11.2—13.6 eV, bottom panel) versus 
-'^^halo- Symbols denote stellar masses of Mstar = 500 Mq (tri- 
angles), 200 Mq (circles) and 25 Mq (crosses). A plateau and a 
small dip at f^^ '^0.1 are due to a weak positive feedback on H2 
formation by an enhanced electron fraction. 



fraction of the gas is neutral but the H2 fraction drops by 
three orders of magnitude from the initial value because 
the LW photons can penetrate further than ionizing ones. 

The evolution after the lifetime of the central star 
{t > 2.2 Myr) is also plotted in Figure El In the 10*^ Mg 
case, if the central star simply fades away without trig- 
gering a supernova explosion (as is assumed in our cal- 
culations), the gas is undisturbed and quickly recombine 
(Fig-Et) and reform molecules (Fig.EJi). The large ion- 
ization fraction within the Hli region can accelerate pro- 
duction of hydrogen molecules through H" formation. 
The shock front continues expanding until it dissipates 
kinetic energy. 

Alternatively, if the central star turns into a super- 
nova, it will generate a blastwave, which runs through 
the ambient gas. A quantity that is of particular im- 
portance in such a case is the final density of the evac- 
uated gas near the center (see 21 for more details). It 
is in fact primarily determined by the host halo mass. 
In the low-mass case (Mhaio = 10^ M©), the central gas 
density drops to '^ 1 cm"^ at i = 2.2 Myr. In the high- 
mass case (Mhaio = 10^ Mq), the central density remains 
much higher at > lO'^cm"'^. The threshold between the 
two cases closely follows that of the escape fraction of 
ionizing photons described in t|3.3l 
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3.2. Photon escape fraction 

A useful output of our simulations is the escape frac- 
tions of ionizing photons and the LW photons from a 
halo. The former is of particular significance in terms 
of early reionization of the Universe (e.g., Yoshida et al. 
2003b; Sokasian et al. 2004; Whalen et al. 2004), while 
the latter quantifies the efficiency of the so-called neg- 
ative feedback on subsequent star formation (Haiman, 
Rees & Loeb 1997; Omukai & Nishi 1999). 

In Figures 0HH1 we plot the escape fractions, averaged 
over the lifetime of the central star, as a function of halo 
mass for various choices of Mstar, w, and Zc. In the case 
of Afstar = 25 Mq (Figure 01 , for instance, the escape 
fraction of ionizing photons /^°" is below 0.5 over the 
range of Afhaio we consider. For Mgtar = 200 Mq, on the 
other hand, f^°^^ is close to unity at Mhaio < 10 Mq but 
drops sharply at Afhaio '^ 2 x 10^ Mq. There appears 
to be a critical mass of the host halo above which /^°" 
is essentially zero. For the escape fraction of the LW 
photons fc^, similar features are found with a slightly 
larger critical mass, because the LW photons are harder 
to shield than ionizing photons. A plateau and a small 
dip at f^^ ^ 0.1 are real; a minor leakage of ionizing 
photons causes a weak positive feedback on II2 formation 
by enhancing electron fraction (Haiman, Rees & Loeb 
1996). 

Figure m indicates that a steeper density profile leads 
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to smaller escape fractions. This is because, for given 
gas mass in a halo, a steeper profile results in higher 
density near the center. The I- front propagation is then 
prevented at early stages due to the higher recombination 
rate. Similarly, an earlier collapse epoch implies a denser, 
compact halo, from which escape fractions are smaller 
(Figure El . 

Strictly speaking, the self-shielding function of the LW 
photons we adopt (Drain & Bertoldi 1996) is derived for 
the stationary gas. In the presence of supersonic flows, 
the photons may doppler-shift out of or into particular 
lines in the LW band. Detailed account for these effects 
requires precise treatment of the line transfer incorporat- 
ing velocity and temperature structures of the gas, which 
is beyond the scope of the present paper. Our results on 
f^^ should therefore be regarded as a minimal estima- 
tion for the true value. 

3.3. Critical mass for photon escape 
Figure \7\ summarizes the degree of radiative feedback 



and photon escape as a function of Aigta 



and 



Low mass halos are ionized promptly after the onset of 
the central star, resulting in high escape fractions of both 
ionizing and the LW photons. On the other hand, the 
Hii regions are heavily confined and the escape fractions 
are essentially zero in high mass halos. For intermediate 
halos, a significant fraction of the LW photons can escape 
from the halos while the ionizing photons are trapped. 
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Figure [3 further indicates that the threshold mass scale 
for the escape of ionizing photons, as a function of Mstar, 
w, and Zc, is well reproduced by the following analytic 
estimation. 

Over the range of gas density gradient we consider 
(1.5 < w < 2.25), the shock velocity is a weak function 
of w and given roughly by 25 ^ 35 km s~^ (see also Shu 
et al. 2002). For the low-mass halos considered in this 
paper, the gravitational infall velocity is 5 ^ 10 km s^^. 
Within the lifetime of the central star iufc, the D-type 
front can then reach the distance 



Id = 20 



t'D 



= 20 f^^^ 

VlOVy V20kms-i 



pc, 



(9) 



where vd is the expansion velocity of a D-type front. 
This is in general much smaller than the virial radius of 
a halo in consideration: 



= 160 



/ Af, 



halo 



1/3 



1 + Zc 

20 



pc. 



(10) 



VIO^Mq, 

applicable to current cosmology at z > 10. It is therefore 
necessary for the I-front to become R-type in order to 
ionize the whole halo. 

As shown in Figure ^ the gas density inside the shock 
front Tg is nearly constant and approximated by the av- 
erage density within Vg: 

3 



-,ii(rg) = 0.39 



1 + Zc 

20 



where ni(rs) is the initial gas density at Tg, normalized by 
the condition that the total gas mass within rvir is equal 
to Ai^haio^e/^M- One can then define the Stromgren 
radius (Stromgren 1939) corresponding to Us as 



r-st = 150 



N^r 



1/3 



1050s- 



Vcm -^ / 



-2/3 



pc. 



(12) 



-1/U)N 



For w > 1.5, Tg depends on rig more weakly (ex Us 

than rst (oc rig ). As far as rst < fg, as is the case 
in the initial stage of expansion, the I-front is of D-type. 
As the shock propagates and Ug decreases, rst eventually 
overtakes rg and the I-front changes into R-type. 

Denoting the radius at which rgt is equal to r^ by roq, 
the necessary condition for the transition into the R-type 



front within the lifetime of the central star is 
namely. 



oq 



where 



M-'H 5.0x10' 



f 'i 



< h, 
(13) 



V270pc 



iVir 



1050s- 



1 + Zc 

20 



3-.^ 



Mr. 



(14) 



Figure 13 shows the above critical mass for complete ion- 
ization, adopting vd = 20 km s~^ in equation ^. In our 

fiducial case with Mgtar = 200 Mg (iVion = 2.6 x lO^o 
s-i, iiifo=2.2 Myr), w = 2, and z^ = 20, it yields 
^CTit — 2.5 X IO^Mq. Indeed it agrees with our sim- 
ulation results over a wide range of parameters. 



If one is to estimate the critical mass for the escape 
of the LW photons separately. Figure H indicates that 
^crit' ~ ^Mjl";" gives a reasonable fit to the simulation 
results for w; ~ 2. In the case of 80 M© < Mgtar S 
500 Mq, 10 < Zc < 30 and w — 2, useful approximations 
for the critical masses are 



M^°"t - 2.5 X 10*= 



Afs, 



200 M, 



M^w ^ 7.5 X 10^ f-^ 
"'* V 200 M( 



3/4 



3/4 



1 



20 



1 + Zc 

20 



-3/2 



M«, 



-3/2 



(15) 



M. 



(16) 

where we have used a power-law fit to the results of 
Schaerer (2002), iVion oc A4\fj, and tufc oc M^°^, over 
the range 80 M© < M^tar < 500 M©. 

If helium is included in the calculation, the gas tem- 
perature will be higher by a factor of '^ 2 for Teff ^ 10^ 
K, because of an increase in the net heating rate (Abel & 
Haehnelt 1999). This will increase the shock velocity and 
consequently lu in equation by ~ 40%. The critical 
mass Afcrit will then become higher by ^ 30% for w = 2. 

3.4. Effective clumping factor of halos 

As an alternative way to quantify the strength of ra- 
diative feedback, we plot in Figure |H1 the ratio Nion/Nn 
as a function of Mgtar, uj, and Zc, where Mon is the total 
number of ionizing photons emitted by the central star 
and A^H is the total number of hydrogen atoms in a halo. 
In order for the ionizing photons to escape from halos 
(/^™ > 1%, circles), Mon/A^H > 10 - 100 is required. 
The threshold values of A'ion/A^H are rather insensitive 
to Afstar but increase with w or Zc. The H2 photodissoci- 
ating photons can escape from halos (/^^ > 1%, crosses 
and circles) for Moh/A'h > 3 ~ 30. 

It is often assumed in the literature that Mon/An needs 
to be higher than the clumping factor {n?) / {n)'^ for the 
gas to be completely ionized. While this assumption is 
likely to be valid for nearly uniform media, it cannot be 
applied to the present case in which the density gradient 
is large. To see this point more clearly, we also display 
in Figure |S1 the effective clumping factor of the halo in 
its initial configuration for 3/2 < w < i: 



Chalo - ("•! )halo/(ni)halo 



(17) 



(3 



w) 



3(2u; - 3) 



1 



X 



(1 



p3-lo^2 ' 



where ( )haio denotes the volume average taken at ri„ < 
'' ^ ^vir and X = 7'in/''vir- The initial condition of the 
present simulations (eq. [7j) gives 



3.93 X 10- 



1 + Zc 

20 



-1 l/to 



(18) 



Over the ranges of parameters studied in this paper, 
C'haio increases with increasing w and decreasing z^ . The 
dependence of Chaio on Zc is simply a product of our 
choice of r-^ (eq. [7])- We have checked that our simu- 
lation results are not sensitive to the specific value of tjh 
(or x), while Chaio is. 

Figure IHl shows that the halos are completely ionized 
and photons can escape from them even for Aion/Ajj < 
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C'haio- The dependences of Chaio on w and Zc do not 
agree with those of the threshold values of Nion/Nn for 
/g™ > 1% either. These are the consequences of the 
enhanced ionization induced by the gas evacuation and 
the subsequent R-type I-front propagation. We stress 
that it is essential to incorporate the dynamical effects 
properly, as in equation (|14(l , to account for ionization of 
the halo with a steep density gradient. 



3.5. Gas evacuation and the evolution of the relic Hll 

regions 

The Hii regions around a massive short-lived star em- 
bedded in a steep density medium never settle into a 
static Stromgren state. The I-front keeps expanding out- 
ward after the central star turns off. The relic Hii re- 
gions then start cooling by recombination and inverse- 
Compton cooling at such high redshifts, and the ionized 
fraction rapidly decreases while the outer part is still ex- 
panding. It is therefore meaningful to investigate how 
much gas is driven out of a halo in the end. In the cos- 
mological context, dark matter halos continue growing, 
so it is likely that the expelled gas eventually falls back 
into the halo by gravitational force. 

Figure El shows the evolution of the baryon fraction 
within Tvir in the cases of Afhaio = 10^ and 10^ Mq in 
our fiducial runs. The halo gas is almost completely 
evacuated by t = 2 Myr for Mhaio = 10^ Mq. For 
Afhaio = 10^ Mq, on the other hand, the gas mass frac- 
tion slightly increases by gas accretion during the first 
few million years. This is simply because the shock ra- 
dius is '^ 70 pc at t = 2 Myr, still inside the virial radius 
(Fig. OJ. Even after the central star fades away, the gas 
will maintain its outward motion until it dissipates ki- 
netic energy. The baryon fraction rapidly decreases at 
~ 5 Myr, when the shock front overtakes the virial ra- 
dius. 

We have further carried out the runs with the same 
set of parameters as those in Figure |51 but omitting the 
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radiation force term in equation ^. While the radia- 
tion force enhances the outward momentum of the gas 
at small radii (Figures Q] and ^ , it helps little to evacu- 
ate the gas from a halo; the baryon fraction within rvir 
is unchanged to within 5% when the radiation force is 
turned off. 

If the central Population III star is so massive that 
they collapse directly to form a blackhole (see, e.g. Fryer, 
Heger & Woosley 2001), it does not disturb the gas me- 
chanically any more. The ultimate fate of the expelled 
gas is then described in simple terms. The cooled gas 
does not provide hydrodynamic pressure in the relic Hii 
region and thus at least some fraction of the gas even- 
tually falls back toward the halo center (i.e., gravita- 
tional potential well) after a while, gaining roughly a 
free-fall velocity. The free-fall time for a halo collaps- 
ing at z = 15 is iff ~ 50 Myr in our adopted ACDM 
cosmology. The host halo itself is growing in mass, but 
the growth timescale is also of the order iff. It is thus 
clear that the subsequent star-formation is possible only 
after this fallback of the gas happens, i.e., gas conden- 
sation at the halo center begins only after a few times 
10^ yr. Note also that the evolution is more complicated 
when the central star explodes as a supernova. In such 
cases, the halo gas may be swept up by the supernova 
blastwave. 

4. COSMOLOGICAL IMPLICATIONS 

4.1. Suppression of star formation 

The internal feedback due to photodissociation of hy- 
drogen molecules by radiation from the very first star 
likely prohibits further star-formation within the same 
halo (Omukai & Nishi 1999). Our results (Figures EHHl 
provide a reasonable estimation for the critical halo mass 
that determines the efficiency of the internal feedback 
in terms of f^^- As shown in §3, early star formation 
via H2 cooling will be self-regulated internally and sup- 
pressed in halos with Mhaio ^ 10^ Mq. For such small 
halos, f^^ is high and the 'one-star-per-halo' assumption 
adopted in the present paper is likely to be adequate. For 
halos with sufficiently larger mass and lower f^^, the 
multiple formation of stars or star clusters could take 
place. Therefore, while the very first stars are likely to 
form in low-mass halos (Abel, Bryan & Norman 2002; 
Bromm, Coppi & Larson 2002), the subsequent star for- 
mation may be suppressed until slightly larger mass sys- 
tems begin to collapse. The latter systems may dominate 
the cosmic star formation rate and metal production at 
high redshifts (e.g. Ricotti et al. 2002). 

Note that our results are not restricted by the under- 
lying one-star-per-halo assumption. We have shown that 
the emitted number of ionizing photons, not the number 
of stars itself, is a key parameter to control the photon 
escape fraction (see eq. d]). Given a weak variation 
of the main sequence lifetime with mass, Mgtar above 
~ 100 Mq can be interpreted as the total mass of stars 
with nearly Eddington luminosities; the case of a single 
500 Mq star is almost equivalent to that of a few 200 M© 
stars. Our results can therefore be rescaled and applied 
to a wider range of stellar numbers. 

4.2. Reionization of the Universe 

The present results on the evolution of Hii regions may 
have profound implications for reionization of the Uni- 



verse. First, high values of photon escape fractions in 
low mass halos imply that such halos can be a major site 
of photon production at least in the early phase of reion- 
ization. In order to surpass recombination at z = 20, 
the emission rate of ionizing photons per unit comov- 
ing volume must be higher than rijon = 3.4 x IO^^Cigm 
s~^ Mpc~'^, where Cigm is the clumping factor of the 
IGM (eq. [26] of Madau, Haardt & Rees 1999). This 
can be achieved if the comoving number of halos host- 
ing a 200 M0 star is greater than 13Cigm Mpc''^. In 
the absence of heavy metals, gas can condense to form 
stars by H2 cooling only in halos above the mass Mh2 — 

5 X lO^/i-i M© (Fuller & Couchman 2000; Yoshida et al. 
2003a). The number of halos just above Mh2, expected 
from the Press & Schechter (1974) mass function, is ~ 30 
Mpc^'^ at z = 20 in the conventional ACDM model with 
{flu, ^A, h, r^B, erg) = (0.3, 0.7, 0.7, 0.05, 0.84). The num- 
ber of low mass halos is thus comparable to that needed 
to ionize the IGM for a moderate value of Cigm- As 
mentioned below, however, the photon production rate 
will be significantly suppressed once secondary feedback 
effects start to operate. 

Secondly, given the sharp decline in photon escape frac- 
tions with increasing halo mass (Figs^0, there remains 
only a narrow range of halo mass in which Population III 
stars can form and contribute to reionization. The pho- 
ton production can therefore be stalled readily as the 
star formation in the low mass halos is suppressed by, 
for instance, the rise in the external UV radiation field 
(e.g., Machacek, Bryan & Abel 2001; Kitayama et al. 
2001; Yoshida et al. 2003a; Oh & Haiman 2003), and 
the mechanical energy input from supernovae (e.g., Ri- 
cotti & Ostriker 2004a). The photon emission may also 
decline as the transition from Population III to Popula- 
tion II stars takes place as a result of metal-enrichment 
(Gen 2003; Mackey, Bromm & Hernquist 2003; Bromm 

6 Loeb 2003). Due to high recombination rate at high 
redshifts, the ionized fraction may even start decreasing 
after a brief period of star-formation in mini-halos. This 
can probably last until a considerable number of stars 
form in larger systems. We will investigate these points 
more quantitatively in our future publication. 

Finally, evacuation of the halo gas by the central star 
may further delay the formation of mini-quasars that are 
suggested to be alternative reionization sources in the 
early universe (e.g. Madau et al. 2004; Ricotti & Os- 
triker 2004b). Such reionization models usually assume 
efficient gas accretion onto central blackholes that are the 
remnants of very massive (Mgtar > 300 M©) Population 
III stars. Our calculations, however, indicate that initial 
gas accretion should be very inefficient for the blackhole 
remnants that formed in small mass (< 10^ M©) halos, 
because the halo gas is effectively evacuated in the first 
place, and, indeed, the gas continues moving outward for 
over 10^ yr after the central star dies. Similar arguments 
hold for the formation of the second generation stars. It 
remains to be seen whether or not early reionization in- 
ferred by the WMAP data is achieved in models with 
these ionizing sources. 

4.3. Supernova feedback 

Mechanical feedback from the first stars is often cited 
as a destructive process in the context of early structure 
formation. If the mass of the central star lies in the range 
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140 M© < Mstar < 260 Mq, it will explode as an ener- 
getic supernova via pair-instability mechanism (Barkat, 
Rakavy, & Sack 1967; Bond, Arnett, & Carr 1984; Fryer, 
Woosley, & Heger 2001; Heger & Woosley 2002). It is 
also suggested that the observed abundance of metal- 
poor halo stars can be accounted for by a "hypernova" 
with progenitor mass 20 Mq < Mgtar < 130 Mq (Umeda 
& Nomoto 2002, 2003). In either case, the energy re- 
leased can be as large as ^ 10^'^ erg, which is much larger 
than the gravitational binding energy of a minihalo with 
mass ~ IO'^Mq. 

At the first sight, this simple argument appears to sup- 
port the notion that high- 2: supernovae are enormously 
destructive. However, the evolution of the cooling super- 
nova remnants (SNRs) crucially depends on the prop- 
erties of the surrounding medium, particularly on the 
central density after the gas is re-distributed by radia- 
tion. Our simulations showed that the final gas density 
at the halo center is primarily determined by the host 
halo mass. For small-mass halos that are almost com- 
pletely ionized by a massive star, the central gas density 
reduces to ~ 1 cm~^. A large fraction of the supernova 
explosion energy will then be converted into the kinetic 
energy of the ambient media. On the other hand, for the 
larger-mass halos, the central gas density exceeds 10^ 
cm~^. Most of the explosion energy will be quickly radi- 
ated away via free- free emission and the inverse-Compton 
scattering of the CMB photons. The free- free emission 
from the hot (T - lO^^^K) SNR leads to soft-Xray emis- 
sion. It has been suggested that positive feedback, in 
terms of molecular hydrogen formation, is possible if an 
early X-ray background builds up (Haiman, Abel & Rees 
2001; Oh'2001; Venkatesan, Giroux & ShuU 2001; see 
however a counter argument by Machacek et al 2003). 
Population III supernova remnants are therefore among 
the most plausible X-ray sources in the early universe. 
We will study in detail the evolution of the SNR and the 
X-ray emission efficiency in a forthcoming paper. 

5. CONCLUSIONS 

We have studied the structure and the evolution of 
early cosmological Hll regions formed around the first 
stars. In particular, we addressed how efficiently the cen- 
tral stars ionize the gas in low-mass halos with Mhaio = 
10^ - 10* M0 collapsing at Zc = 10 - 30. We showed 
that a single massive star can ionize all the halo gas in a 
few million years if the host halo mass is smaller than a 



few million solar masses. The final ionization fraction of 
the halo gas depends sensitively on the initial gas density 
profile as well as the host halo mass. While a few pre- 
vious cosmological simulations suggest an approximately 
isothermal density profile [p oc r~^) around the first star- 
forming regions, the exact slope may vary among the pri- 
mordial gas clouds. We therefore carried out a number 
of simulations for plausible cases and explored a large 
parameter space. For a similar set of parameters, our 
results are in good agreement with those of Whalen et 
al. (2004), who studied the case of Mhaio = 7 x 10^ Mq 
and Zc ~ 20. 

The formation of the Hii region is characterized by 
initial slow expansion of a weak D-type ionization front 
near the center, followed by rapid propagation of an R- 
type front throughout the outer gas envelope. We find 
that the transition between the two front types is in- 
deed a critical condition for the complete ionization of 
halos of cosmological interest. This accounts for the fact 
that the photon escape fraction has a sharp transition 
from ^ 100% to ^ 0% at a certain critical mass scale of 
10^ ^ 10^ Mq. It is also responsible for the fact that the 
whole halo can be ionized under smaller values of photon- 
to-baryon ratio than the effective clumping factor. The 
radiation force can contribute to enhancing the outward 
motion within ^ 30 pc at the initial stage of expansion. 

Based on the obtained numerical results, we developed 
an analytic model to predict the escape fractions of ioniz- 
ing photons and dissociating photons. We find that there 
is a narrow range of halo mass in which Population III 
stars can form and contribute to ionizing the Universe, 
provided that low mass halos are in infall phases and host 
a small number of massive stars. The analytic model and 
the numerical results presented in this paper provide use- 
ful ingredients for the studies of early reionization, the 
evolution of high-rcdshift supernova remnants, and the 
formation of the second generation stars. 
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APPENDIX 
IMPLEMENTATION OF RADIATIVE TRANSFER 

Above the Lyman limit of hydrogen [hv > 13.6 eV), we solve the following radiative transfer equation taking account 
of both absorption and emission: 

-T- = -nm<yvlu + jy, (Al) 

as 

where I^ is the specific intensity at frequency v, s is the distance along a light ray, Ui, is the photoionization cross 
section of hydrogen (Osterbrock 1989), and ji, is the local emission coefficient for ionizing radiation. In the case of 
kT <ti hi>h as in the present simulations, ji, is given by (Tajiri & Umemura 1998; Susa & Umemura 2000) 



Jy 



hu a I Tie null 
47rdV 








5v, 



(A2) 



where ui^ is the Lyman limit frequency, ai is the recombination rate to the ground state of hydrogen (Osterbrock 1989), 
and Sv = k T /h is the thermal width of the recombination line. Recombinations to the excited states are excluded 
in eq uation ljA2|) because they are unable to produce ionizing photons above the Lyman limit. In practice, equation 
ljAl|l is solved separately for vi, < v < ui^ + 5v and v > i/]^ + 5v. For the latter, photoionization is regarded as pure 
absorption and the solution reduces to 



Iy{s) = 1^(0) cxp (-fTixA^Hi) , {y>vi, + 8v) 



(A3) 



where A/hi = /g TiYiids' is the column density of neutral hydrogen along the photon path. 

The coefficient (per unit time) and the heating rate (energy per unit time per unit volume) for photoionization of 
hydrogen atoms (H + hv — > e^ + H"*") are then given respectively by 
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(A4) 
(A5) 



where the integral with respect to fi is taken over all solid angles. Assuming that the momentum of a photon is 
transfered to the gas upon each ionization, the radiation force per unit mass in the direction specified by a unit vector 
n is expressed as 

/rad(n^) — dfl dv I^ayCOsO (A6) 

"IpC J J^^ 
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where 9 is the angle between n and the hght ray, m-p is the proton mass, and c is the speed of hght. In our spherically 
symmetric simulations, n lies in the radial direction. Photon momentum carried away by reprocessed radiation is 
properly taken into account by using the solutions of multi-frequency radiative transfer mentioned above. 

As far as only the radiation from a central star is concerned, the intensity integrated over solid angles is well 
approximated at v > vi^ -\- 5v by 

7- JO [t a^n ^1^ exp(-cr^A/'Hi) 

Ii,dil ^ I,ycos0dil— ^ , (Ay) 

J Airr'^ 

where L^ is the stellar luminosity, and r is the radial distance from the center, taken to be sufficiently larger than the 
physical size of the star. On the other hand, the diffuse radiation just above the Lyman limit, j^l < i^ < J^l + ^v, comes 
from all the directions. For the diffuse component, we directly solve equation (|Aip by means of an impact parameter 
method (Mihalas & Weibel-Mihalas 1984). At each radial point, angular integrations in equations (|A4l) - (|A6p are done 
over at least 20 bins in 6* = — tt. This is achieved by handling 400-1,000 impact parameters for light rays. 

In the LW (11.2-13.6 eV) band, we use the self-shielding function of Draine & Bertoldi (1996) to compute the 
coefficient (per unit time) for photodissociation of II2 molecules (H2 + /i;^ ^ H + H) : 

fcH2 = O.lS^pumpX/shiold, (A8) 

where 

epu„,p~3xl0-i°s-i, (A9) 

4 X lO^^"' erg cm^'^ 4 x 10~^^ erg cm 

k/H2 



X= ..,,n^l. -3 - ..,.n-14 -3 ' (^10) 



1, 



1014 cm-2 



/shield = min 

and 7\/h2 is the column density of H2 
We also compute the coefficients f 
(HJ + hv ^ 11+ + H) from the radiation below vi^ as 



A/h 



(All) 



We also compute the coefficients for photodetachment of H (H -\- hv ^r e +11) and photodissociation of Hr'" 



2 



where the subscript i stands for H^ or Hj and a^^i, denotes the cross section. For photodetachment of H", we use a 
fitting formula by Tegmark et al. (1997) to the cross section of Wishart (1979) with the low-energy cut-off hv^i- — 0.74 
eV. For photodissociation of H^, we use the cross section from Table 2 of Stancil (1994) down to hv^+ = 0.062 eV. 

Our results are insensitive to specific values of the low-energy cut-off, because black-body spectra of stellar sources 
diminish at low energies. The radiation above vi^ is also irrelevant since it will be heavily absorbed by hydrogen atoms 
whenever H~ and H^ processes become relevant. 

CODE TEST 

We have tested and verified the accuracy of our Lagrangian code by a variety of test problems, such as the shock-tube 
problem, the Sedov explosion problem (Sedov 1959), and a comparison with the self-similar solution for the adiabatic 
accretion of coUisional gas (Bertschinger 1985). We here describe the results of a "Stromgren sphere" test, which is 
the most relevant to the present paper. 

We simulate the propagation of an I-front around a central source with iVjon = 3.3 x 10^° s^^ and Tcff — 10* K in 
a uniform medium with rifj = 1 cm~^. The Stromgren radius is at rgt — 220 pc for T — 10^ K. For the sake of direct 
comparison with analytic solutions, we exclude the radiation force and gravitational force terms from the momentum 
equation (eq. |2]). Figure lBlOl shows the radial profiles of hydrogen density, temperature, HI and H2 fractions, and the 
gas velocity. In contrast to the cases with steep density gradient {w > 3/2), the I-front is initially R-type at r < rgt 
and then changes into D-type at r > rgt. The I-front is resolved with more than 10 gas shells. As in Figure ^ there 
temporally appears a thin shell of II2 just in front of the Hil region. 

In Figure IBIII we also plot the evolution of the I-front radius. For definiteness, we define it as the radius at which 
HI fraction is equal to 0.1. Also plotted for reference are the analytic solutions of I-front propagation in stationary 
and dynamical media (eqs [12-9] and [12-20] of Spitzer 1978) for the photoheated temperature of T = 10^ K. They are 
expected to describe the evolution of R-type and D-type fronts, respectively. Our simulation results in fact reproduce 
both of them within 5% accuracy at i > 10'^ yr. This further ensures the capability of the code for the simulations in 
the present paper. 
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Fig. bio. — Propagation of an I-front in a uniform medium around a source with A^ion = 3.3 X 10^ s~ and T^ff = 10 K. Radial profiles 
are shown at t = (dotted lines), 10^ yr (dashed lines), and 10^ yr (solid lines) for (a) hydrogen density, (b) temperature, (c) HI (thin 
lines) and H2 (thick lines) fractions, and (f) radial velocity. The Stromgren radius is at rst = 220 pc. 
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Fig. Bll. — Evolution of an I-front in the simulation shown in Figure lb 101 ('circles') . Lines denote the analytic solutions for the I-front 
propagation in stationary (solid) and dynamical (dashed) media described in the text. 



